# 加载相关包
library(openxlsx)
library(geepack)
# 加载自定义函数
source("./func/gee_run.r")

############# 需求 ############# 
# 1. 多因素删除缺失值：分别剔除因变量和自变量中有缺失值的行，然后做多因素 gee 分析
# 2. 单因素虚拟填补缺失值：使用虚拟值填充缺失值，然后做单因素 gee 分析
# 3. 多因素 mice 填补缺失值：使用 mice 填补因变量和协变量的缺失值，然后做多因素 gee 分析
############# 需求 ############# 

###### 读取数据，数据预处理
data_test = read.xlsx("./input/data918.xlsx", 3)
###
data_tmp = data_test[data_test$month != 6, ]
data_tmp = data_tmp[order(data_tmp$id, data_tmp$month), ]
extract_id = unique(data_tmp[duplicated(data_tmp[, 1]), ]$id)
data_tmp = data_tmp[data_tmp$id %in% extract_id, ]
###

data_anal = data_tmp[, c(1, 9:13, 14:25, 27, 29:39)]
varsx = colnames(data_test)[c(14:25, 27)]  # 自变量
varsy = colnames(data_test)[c(9:13)]       # 因变量
varsc = colnames(data_test)[c(29:39)]      # 协变量

# 确保变量类型一致：因变量为数值，自变量、协变量为因子
data_anal$id = as.character(data_anal$id)
data_anal[varsy] = lapply(data_anal[varsy], as.numeric)
data_anal[varsx] = lapply(data_anal[varsx], as.factor)
data_anal[varsc] = lapply(data_anal[varsc], as.factor)
# str(data_anal)
######


# 1. 多因素删除缺失值
gee_run(data_anal,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc[-5],
        file = "./output/0918_delete_10cv_del_birthweight_c.xlsx")
gee_run(data_anal,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc[-2],
        file = "./output/0918_delete_10cv_del_week_c.xlsx")



# 2. 单因素虚拟填补缺失值
## 2.1 预处理，将缺失值定义为该列的最大值 + 1
tmp_data = data_anal
tmp_data[varsx] = lapply(tmp_data[varsx], as.numeric)  # 先转换为数值，才能取得到最大值
for (k in varsx) {
  tmp_data[, k] = ifelse(is.na(tmp_data[, k]), max(tmp_data[, k], na.rm = T), tmp_data[, k])
}
# tmp_data[is.na(tmp_data)] = 99  # 或者快捷设置为一个特定值
tmp_data[varsx] = lapply(tmp_data[varsx], as.factor)   # 再转换回因子
## 2.2 gee（单因素不要传入 varsc 参数）
gee_run(tmp_data,
        varsx = varsx,
        varsy = varsy,
        file = "./output/0918_virtual.xlsx")


tmp_data = data_anal
tmp_data[varsx] = lapply(tmp_data[varsx], as.numeric)  # 先转换为数值
tmp_data[varsc] = lapply(tmp_data[varsc], as.numeric)
tmp_data[is.na(tmp_data)] = 9
tmp_data[varsx] = lapply(tmp_data[varsx], as.factor)   # 再转换回因子
tmp_data[varsc] = lapply(tmp_data[varsc], as.factor)
## 2.2 gee（单因素不要传入 varsc 参数）
gee_run(tmp_data,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc,
        file = "./output/0922_virtual.xlsx")


# 3. 多因素 mice 填补缺失值
## 3.1 先跑一次 mice
tmp_data = unique(data_anal)
miss = md.pattern(tmp_data)                # 查看缺失数据的模式
imp = mice(tmp_data, m = 5, seed = 1)      # 填补
mice_data = complete(imp, action = 1)      # 懒惰的人选择第一个
write.xlsx(mice_data, "./output/output_mice_data.xlsx", overwrite = T)   # 写入 mice 后的数据
## 3.2 gee
gee_run(mice_data,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc[-5],
        file = "./output/0918_mice_10cv_del_birthweight_c.xlsx")
gee_run(mice_data,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc[-2],
        file = "./output/0918_mice_10cv_del_week_c.xlsx")
